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Combining a semi-classical analysis with exact diagonalizations, we show that the ground state of the SU(3) 
Heisenberg model on the square lattice develops three-sublattice long-range order. This surprising pattern for 
a bipartite lattice with only nearest-neighbor interactions is shown to be the consequence of a subtle quantum 
order-by-disorder mechanism. By contrast, thermal fluctuations favor two-sublattice configurations via entropic 
selection. These results are shown to extend to the cubic lattice, and experimental implications for the Mott- 
insulating states of three-flavor fermionic atoms in optical lattices are discussed. 
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Mott transitions, i.e. metal-insulator transitions driven by 
correlations, and the nature of the associated Mott insulating 
phases represent one of the central themes of contemporary 
condensed matter physics [1], and more recently also of the 
field of ultracold atomic gases [2]. Theoretically, the canoni- 
cal case of two-flavor fermions on hypercubic lattices is thor- 
oughly understood. For strong interactions an antiferromag- 
netically ordered two-sublattice Neel state is realized. Ongo- 
ing experimental efforts using ultracold fermionic gases are 
focused on reaching this state coming from higher tempera- 
tures [3]. 

In an exciting parallel development, recent experimental 
advances using multi-flavor atomic gases [4, 5] have paved the 
way to the investigation of Mott insulating states with more 
than two flavors in optical lattices [6, 7]. While it is intu- 
itively clear that Mott insulating states will exist at particular 
commensurate fillings - as suggested by atomic limit consid- 
erations and single site DMFT simulations [8, 9] -, the nature 
and the spatial structure of multi-flavor Mott insulating states 
are in general not well understood. For instance, on the square 
lattice geometry, many different proposals for insulating states 
have been put forward, ranging from SU(iV) symmetry break- 
ing "magnetic" states to dimerized or plaquette states, chiral 
spin liquids and staggered flux phases [6, 10-13]. 

In this Letter, we present a strong case in favor of a three- 
sublattice long-range ordered ground state for the Mott insu- 
lating state of three-flavor (N = 3) fermions with one particle 
per site ( 1/3— filling) on the square lattice. This is based on an- 
alytical and numerical investigations of the strong coupling 
limit U ^ t of the SU(3) symmetric Hubbard model defined 
by the Hamiltonian 

H=-tJ2 {c[ a c ja +h.c.) +U "i.e.ni./J ■ (1) 

(i,j),a i,a<P 

Here c\ a and c i a create and annihilate a fermion at site i with 

flavor a, respectively, and n^ a = cj Q c ; a . To second order in 
t/U, the low-energy physics is captured by the SU(3) antifer- 




FIG. 1: (Color online) Sketch of the (a) three-sublattice and (b) two- 
sublattice phases of the SU(3) antiferromagnetic Heisenberg model. 



romagnetic Heisenberg model with coupling J = 2t 2 /U: 

where Vij is a transposition operator that exchanges SU(3) 
spins on site i and j: Vi^otifij) = The spins on a site 

form the 3-dimensional fundamental irreducible representa- 
tion of the SU(3) algebra. In the following, the basis states 
will be denoted by \A), \B) and |C). Note that this model can 
also be seen as a special high-symmetry point of the SU(2) 
spin-1 bilinear-biquadratic exchange Hamiltonian when bi- 
linear and biquadratic couplings are equal. 

Let us start by briefly reviewing what is known about this 
model in various geometries. In one dimension, the model 
has a Bethe Ansatz solution [14]. It has gapless excitations at 
q = and q = ±27r/3 [15], and the correlations decay al- 
gebraically with period 3. In higher dimensions, much less is 
known, and most of it relies on the pioneering work of Papan- 
icolaou [16] who has investigated this question in the context 
of spin-1 models with the help of a variational approach, us- 
ing a site-factorized wavefunction of the form 

N A 

l*> = II ( d A,i\ A )i + d B, l \B) l + dcAC)i) , (3) 
1=1 

where Na is the number of sites. Grouping the variational 
parameters into (complex) vectors d; = (d>A,i, dB,i,dc,i) and 
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imposing the normalization • = 1, the problem reduces 
to the minimization of 



E v 



(*\H\*) 



J^ldi-d.-l 



(4) 



Since J > 0, the energy of nearest-neighbor bonds is minimal 
when the d vectors, hence the wavefunctions, are orthogo- 
nal. The nature of the variational ground state then depends 
on the connectivity of the lattice. For a triangle, this condition 
enforces three mutually orthogonal d vectors. As a conse- 
quence, for the triangular lattice, the energy is minimized up 
to global SU(3) rotations by a single wavefunction constructed 
by choosing on each of the three sublattices the wavefunctions 
\A), \B), and |C), respectively. The three-sublattice long- 
range order embodied by this wavefunction has recently been 
shown to be stable against quantum fluctuations [17]. 

By contrast, the square lattice does not provide enough 
constraints to uniquely select a set of mutually orthogonal 
d vectors. Consider for instance a Neel state with \A) and 
\B) on the two sublattices. Any state obtained by replacing 
\B) by |C) on an arbitrary number of sites is also a ground 
state. This variational approach thus leads to a highly de- 
generate ground-state manifold. This situation is reminiscent 
of frustrated SU(2) antiferromagnetism, where the competi- 
tion between exchange paths often leads to an infinite num- 
ber of classical ground states: in that case, quantum or ther- 
mal fluctuations often restore long-range order by a selection 
mechanism that favors collinear or planar configurations and 
is known as 'order-by-disorder' [18]. For the SU(3) model, 
zero-point quantum fluctuations can be calculated with the 
help of the flavor-wave theory, an extension of the SU(2) spin- 
wave theory to the SU(3) case [16, 19]. This approach starts 
from the representation of the model in terms of three-flavor 
Schwinger bosons: 



H,u&{A,B,C} 
,t 



(5) 



with the constraint al i a v ,- L = 1. In order to treat quan- 
tum fluctuations around a variational solution defined by dj, 
one first performs a local SU(3) rotation of the Schwinger 
bosons by choosing two vectors ej and fj which, together 
with d 3 , define a local orthogonal basis, in terms of which 
the rotated bosons are given as a\ ■ — ^ d^ja^j, d B ■ = 

Eju e M,J a L'' md a c,j = E^/w^j. A 'semiclassical' 
1 /M expansion is then generated by replacing the constraint 

by Si/ a t i a v i = M and the Schwinger bosons along the lo- 
cal direction of the variational solution by 



a\,nO>A,i \/ M -an_,d 



st 7, ~J 7, 



1 'B,i u 'B,i a C,i a C,i 



(6) 



The bosons a,B,i and dc,i play the role of Holstein-Primakoff 
bosons. Expanding in powers of 1/M leads to 

Vij = M (dj ■ Aj + d, ■ Aj\ (dj ■ A i + d 4 ■ A]) -M (7) 



with Aj = eid Bi + fid Ci . The resulting Hamiltonian is 
quadratic and can be diagonalized by a Bogoliubov transfor- 
mation. The zero-point energy is half the sum of the eigenfre- 
quencies. 

Now, simple considerations lead to two natural candidates 
for possible orderings: (i) the two-sublattice state with order- 
ing wavevector Q2 = (7r,7r), as suggested by the bipartite 
nature of the lattice; (ii) a three-sublattice state with ordering 
wavevector = (2ir/3, ±27r/3), as suggested by the soft- 
ening of the excitation spectrum in one-dimensional chains. 
Both states feature diagonal stripes of equal spins, with an al- 
ternation of A and B stripes in case (i) and a succession of A, 
B and C stripes in case (ii) (see Fig. 1). For such states, the 
excitation spectrum has the periodicity of the square lattice, 
leading to 

2 / 1 

n = -2M JN A + M^2^2cj v (k) ( a£(kK(k) + - 

u=l k ^ 

(8) 

The spectrum consists of two branches because there are two 
states orthogonal to the local variational state. For the three- 
sublattice state, the branches are degenerate with dispersion 



w li2 (k) = 23yf\ - | 7 k| 2 , (9) 
while for the two-sublattice case they are given by 



Wl (k) = 2jA/4-(7k + 7kr , wa(k)=0, 



(10) 



with 7k = (e 4fcx + e lky )/2. The corresponding spectra are 
shown in Fig. 2(a) and (c), respectively. 

In both cases, there is an infinite number of zero modes. For 
the two-sublattice case, the presence of an entirely soft branch 
is natural: if the two-sublattice variational state corresponds, 
say, to A and B spins, each spin can tilt toward C at no en- 
ergy cost. What is more surprising is the line of zero modes in 
the three-sublattice case, since one would naively expect zero 
modes only at k = 0, ±Q^. This suggests that there must be 
a family of classical ground states with wave vectors (q, ±g). 
Indeed, after choosing the vectors d; and d; + i of two consec- 
utive stripes, we may write the next one as d; + 2 = cos 9di + 
sin^df x d; + i. If we continue this construction using al- 
ways the same 6, we end up with a helical state of the form 
\ipi) oc y/- cosq\u)+smq(x l +y l )\v)+cosq(x t +y i )\w) with 
cos9 = — 1 — 2 cos q and ir/2 < q < ir, where \u), \v), \w) 
depend on the initial choice of d; and d; + i. These helical 
states extrapolate between the two-sublattice case (q = n/2 
and q = ir) and the three-sublattice case (q = 2tt/3). Their 
excitations have the periodicity of the lattice, and the disper- 
sion is given by the non-negative solutions of the equation 
(extended to a cubic lattice of arbitrary dimension D for later 
reference) 



D 2 J 2 [2(1 - 7k7k ) + (2 - 72 - 7 2) cos 2 9] uj 2 + 
0, (11) 



D 4 J 4 sin 4 9(1 - 7k7k ) 2 
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FIG. 2: (Color online) Flavor-wave dispersions of (a) the three-sublattice state, (b) the helical state with cos 9 = 0.4, and (c) the two-sublattice 
state, (d) Zero-point energy per site, ezp = (1/JVa) X)k v ^"(10/2, as a function of the wavevector (q, ±g) of the helix. The minimum is 
located at q — 2tt/3 in both two and three dimensions. 



where 7k = (e * + e lky + . . .)/D. An example is shown in 
Fig. 2(b). All helical states give rise to a line of zero modes. 

We see therefore that there is in fact an infinite number of 
helical ground states that are a priori good candidates for be- 
ing the quantum ground state. We have calculated and com- 
pared the zero-point energy of these states (Fig. 2), and found 
that it is minimal for q = 2ir/3. We have also compared the 
zero-point energy of the three-sublattice state with that of ran- 
dom ground states on finite clusters, as well as with that of 
stripe states with lower periodicity, with the conclusion that 
it is always lower. So, within the flavor-wave theory, the sta- 
bilization of the three-sublattice state appears to be quite ro- 
bust. This leads to the first important conclusion of this pa- 
per: on the square lattice, quantum fluctuations stabilize three- 
sublattice long-range order. 

This conclusion is quite surprising from the point of view of 
order-by-disorder. Indeed, the two-sublattice state has by far 
the largest number of zero modes, and according to common 
wisdom it should be selected. However, this need not be the 
case for quantum fluctuations: if the non-zero modes have suf- 
ficiently large energy, they may compensate for the vanishing 
contribution of the zero modes. This is what happens here for 
the two-sublattice structure, whose upper branch is larger than 
twice the degenerate three-sublattice branch for all wavevec- 
tors. By contrast, for thermal fluctuations, the stabilization 
of the configuration with the largest number of zero modes is 
systematic, since the low-temperature free energy reads [20] 

F = E - —^-TkiT — Nm NzM TlnT , (12) 
4 2 

where Nm is the total number of modes and Nzm is the num- 
ber of zero modes. The classical spectrum A is found by di- 
agonalizing the quadratic form obtained by replacing creation 
and annihilation operators by complex numbers. For the heli- 
cal states, its four branches are given by 

\ 2 —DJ [2 ± cos 0(7k + 7k)] X+D 2 J 2 sin 2 6(1- 7k 7k) = . 

(13) 

Similarly to the SU(2) case, the classical and quantum spectra 
are different but related, and in dimension D > 2, the propor- 
tion of zero modes is the same. This analysis predicts there- 



fore that thermal fluctuations stabilize the two-sublattice state, 
in agreement with Classical Monte Carlo simulations [21]. 

Next, we compare these predictions with exact diagonaliza- 
tions (ED) of the Hamiltonian of Eq. (2) on finite clusters. The 
energy per site for square samples of up to 20 sites is shown 
in Fig. 3(a). It is significantly smaller for the samples whose 
number of sites is a multiple of 3 (9 and 18), providing ev- 
idence in favor of a three-sublattice symmetry breaking. To 
check if the continuous SU(3) symmetry is also broken, we 
have plotted in Fig. 3(b) the energy levels as a function of the 
quadratic Casimir operator C\ of SU(3), keeping track of the 
irreducible representations of the space group symmetry. In 
the case of the continuous symmetry breaking one expects the 
low-energy part of the spectrum to align linearly as a func- 
tion of C2, giving rise to a tower of states [22, 23]. This is 
clearly the case in Fig. 3(b), as highlighted by the dashed line, 
This tower of states can be thought of as a combination of two 
towers corresponding to the two possible propagation vectors 
(Q3 ), which results into the finite-size splitting of some lev- 
els (e.g. TA1 and TB2), as well as the increased degener- 
acy of some irreps (e.g. W). Note also that the tower is not 
as well separated from other states as in other systems [22], 
a consequence of the order-by-disorder selection mechanism 
that leads to low-lying excitations associated to other mean- 
field solutions. The structure of the energy spectrum further 
indicates that the state with an equal population of the SU(3) 
basis states is stable with respect to the occurrence of spon- 
taneous population imbalance or phase separation. Finally, 
an inspection of the real-space correlation functions of the 18- 
site sample allows for a rough estimate of the ordered moment 
of about 60-70% of the saturation value. Long-wavelength 
fluctuations on larger systems might further reduce this value, 
but with such a large value on 18 sites, we expect the or- 
der to survive in the thermodynamic limit. Note that, due 
to zero modes, estimating the ordered moment within flavor- 
wave theory would require pushing the expansion beyond lin- 
ear order. So altogether, ED provide very clear evidence in 
favor of the three-sublattice flavor-wave state. 

Let us now briefly discuss the experimental implications of 
these results. Reaching sufficiently low temperatures is cur- 
rently a major challenge in ultracold atomic systems : inter- 
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Finally, the detection of the three-sublattice order might be 
attempted using noise correlations [25], since the structure 
factor is expected to have a peak at . Alternatively, a recent 
report of single atom resolution experiments [26] suggests that 
direct imaging might be possible provided that some contrast 
can be achieved between different atomic species. 

We thank M. Hermele and V. Gurarie for interesting discus- 
sions on this subject. This work was supported by the Swiss 
National Fund, by MaNEP and by the Hungarian OTKA Grant 
No. K73455. F.M. and K.P. are thankful for the hospitality of 
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FIG. 3: (Color online) (a) Energy per site from exact diagonaliza- 
tion of various finite size square clusters compared to the flavor-wave 
(FW) result, (b) Tower of states for 18 sites. T denotes points in 
the center of the Brillouin zone, with Al and A2 one-dimensional 
and El two-dimensional irreps, while W is a four dimensional irrep 
with wave vectors (±27r/3, ±27r/3). The tower contains the irre- 
ducible representations expected for a striped three-sublattice flavor- 
wave state with two possible orientations. 

estingly, the exchange integral of the SU(N) case is equal to 
2t 2 /U, independently of N. It is thus realistic to expect that 
the exchange scale can be reached for SU(N) fermions as soon 
as it is reached for SU(2) ones. In that respect, it will be im- 
portant in experiments to carefully choose the optimal cou- 
pling strength U/t, which should be large enough to put the 
system into the Mott insulating phase described by the SU(3) 
Heisenberg model, but not too large to lead to accessible val- 
ues of the energy scale set by the exchange integral. 

Next, we note that the conclusions regarding the selection 
by quantum or thermal fluctuations for SU(3) fermions are ac- 
tually valid in any dimension D > 2, as can be checked easily 
fromEqs. (11,13). So the present results allow us to make pre- 
dictions both for the square and cubic lattices, and the compe- 
tition between quantum and thermal order-by-disorder should 
lead to a rather rich physics. In both 2D and 3D, the system 
should first develop two-sublattice ordering tendencies as it 
is cooled below the exchange scale. In two dimensions, we 
expect the system to undergo a finite temperature transition 
at lower temperatures into a directionally ordered state (se- 
lection between the two independent spiral propagation 
vectors), and to reach a three-sublattice ordered state at zero 
temperature. In three dimensions however, a finite tempera- 
ture transition into a two-sublattice ordered state is in prin- 
ciple possible, leading to two possible scenarios: upon low- 
ering the temperature, the system might first undergo a tran- 
sition into a two-sublattice ordered state, which is followed 
by a second transition into the three-sublattice ordered state 
selected by quantum fluctuations. Alternatively, it could un- 
dergo a direct first-order transition from the paramagnetic into 
the three-sublattice ordered state. High-temperature series ex- 
pansion of the SU(N) case on the 3D cubic lattice seems to 
favor the second possibility [24]. 
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